/*
Use 8000 results from the server to pick the calibration
*/
/*
This program chooses one of 8000 calibrations on the server
*/
clear
tempfile temp
tempfile temp2

set more off
set scheme sonia_new
local exp all_experience_l1

************************************************************************************************************************************************************************
************* Import Results From Data (saved in save_calibrationcont`exp'.xlsx, made in calibration_static_continuous.do in the analysis folder) **********************
************************************************************************************************************************************************************************
* Need exit rates, entry rates, distribution of agents (for weights), listing distribution(?) take from the appendix .77 (num_bs on experience)

*set up data experience 1-100, for each 9 states
set obs 900
gen state = floor((_n-1)/100)+1
gen experience = mod(_n,100)
replace experience=100 if experience==0 
tempfile temp
save `temp'


import excel data_for_model/save_calibrationcont`exp'.xlsx, sheet("states") clear // entry rates for each 9 states, what a horrible name.. 
gen state = (floor(A/10)-1)*3+mod(A,10)
rename B entry_data
drop A
merge 1:m state using `temp'
drop _merge
tempfile temp
save `temp'


import excel data_for_model/save_calibrationcont`exp'.xlsx, sheet("states_exit") clear // exit rates for each 9 states, what a horrible name.. 
xpose, clear
rename v1 state_11
rename v2 state_21
rename v3 state_32
rename v4 state_33
drop if _n == 1
gen experience = _n
reshape long state_, i(experience) j(A)
gen state = (floor(A/10)-1)*3+mod(A,10)
rename state_ exit_rate_data
drop A
sort state experience
drop if experience>100

merge 1:1 state experience using `temp'
drop _merge

gen b_data = .77
tempfile temp
save `temp'


foreach state in bust medium boom {
	import excel data_for_model/save_calibrationcont`exp'.xlsx, sheet("`state'") clear // sale probability for each of the three states
	xpose, clear
	keep v1 v4 v8
	rename v1 experience
	replace experience = experience+1
	drop if experience>100
	rename v4 sale_prob_data
	rename v8 learning_data
	expand 3
	bysort experience: gen state = _n
	if "`state'"=="bust" {
		replace state = (state-1)*3+1
	}
	if "`state'"=="medium" {
		replace state = (state-1)*3+2
	}
	if "`state'"=="boom" {
		replace state = (state-1)*3+3
	}
	
	merge 1:1 state experience using `temp'
	drop _merge

	tempfile temp
	save `temp'
}

import excel data_for_model/save_calibrationcont`exp'.xlsx, sheet("distribution_state") clear // distribution of agents across states
gen state = (floor(A/10)-1)*3+mod(A,10)
drop A
rename B experience
replace experience = experience+1
rename C dist_data
bysort state: egen total_agents_data = total(dist_data) if experience<=50
gen density_data = dist_data/total_agents_data if experience<=50
merge 1:1 state experience using `temp'
drop _merge
tempfile temp
save `temp'


tempfile data_moments
save `data_moments'
save data_moments, replace 
**********************************************************************************************************************************
*** Import Results From Model (compiled in compile_results.m that combines all the .mat files for every grid point ) ************
**********************************************************************************************************************************


import delimited using calibration_results_data.csv, clear

rename v1 experience
rename v2 state
rename v3 listings_model
rename v4 buyers_model
rename v5 sale_prob_model
rename v6 profit_model
rename v7 exit_rate_model // exit rates
rename v8 dist_model
rename v9 entry_model
rename v10 sellV_model
rename v11 BuyPr_model
rename v12 mu
rename v13 sigma
rename v14 c_e
rename v15 delt
rename v16 b_model // listing accumulation
rename v17 failed // listing accumulation
drop if failed==1
gen num_bs_model = listings_model +buyers_model

* check there is no same calibration twice
drop if missing(mu) | missing(sigma) | missing(c_e) | missing(delt) | missing(state) | missing(experience)
duplicates drop
bysort  mu sigma c_e delt state experience: gen i = _n
keep if i==1

isid mu sigma c_e delt state experience 
egen modelid = group(mu sigma delt c_e)
sort modelid state experience

* This is stationary distribution of states, we should weight the results according to this OR 
* according to states we actually observe in the data
gen q = 0
replace q = 0.2209 if state==1
replace q = 0.052 if state==2
replace q = 0.065 if state==3
replace q = 0.075 if state==4
replace q = 0.1864 if state==5
replace q = 0.0621 if state==6
replace q = 0.0424 if state==7
replace q = 0.0847 if state==8
replace q = 0.2119 if state==9

gen bust=1 if state==1 | state==4 | state==7
gen medium=1 if state==2 | state==5 | state==8
gen boom=1 if state==3 | state==6 | state==9 


* entry rates
bysort modelid state: egen num_agents_model = total(dist_model)
gen entry_rate_model = entry_model/num_agents_model

* density
gen  density_model = dist_model/num_agents_model // will be useful for model weights

merge m:1 state experience using data_moments.dta
drop if experience>50
drop _merge

**********************************************************************************************************************************
*********************** Compute moments and total optimization value for each parameter combination ******************************
**********************************************************************************************************************************

* compute moments
gen mom_exit = abs((exit_rate_model-exit_rate_data)/exit_rate_data)
gen mom_entry = abs((entry_rate_model-entry_data)/entry_data)
gen mom_list = abs((b_model-b_data)/b_data)
gen mom_dist = abs((density_model-density_data)/density_data)

* compute weights
gen weight = max(0,q*density_model) 

* weighted moments
bysort mu sigma c_e delt: egen wt_mom_dist = wtmean(mom_dist^2), weight(weight)

* minimization functions
gen total_moment = mom_dist^2 + mom_exit^2 + mom_entry^2
bysort modelid: egen val_fun = wtmean(total_moment), weight(weight)
sort val_fun state experience
sum val_fun // this is where I look up the parameters that match the minimum va_fun

local minimum = `r(min)'
sum mu if (val_fun-`minimum')<0.00001
sum sigma if (val_fun-`minimum')<0.00001
sum c_e if (val_fun-`minimum')<0.00001
sum delt if (val_fun-`minimum')<0.00001
